You don't actually need to make a new PR each time you make changes. If you make further commits on the same branch and push to your GitHub account, it automatically adds those commits to the previous PR.

I think you are correct. Hmm, I guess the test code that I wrote does not check this :)

Just checking: this PR replaces the previous one, right?

Don't know about the SNF. One easy is to repeatedly do HNF of the matrix and its transpose.

Some other references for efficient sparse integer matrix algorithms:
"Algorithms for Large Integer Matrix Problems", by Mark GiesbrechtMichael Jr. JacobsonArne Storjohann:
https://link.springer.com/chapter/10.1007/3-540-45624-4_31
and "A note on the hermite basis computation of large integer matrices" by Vollmer : https://dl.acm.org/doi/10.1145/860854.860905

There is also https://cs.uwaterloo.ca/~astorjoh/iml.html (in case you did not know).

Sorry, I am not an expert on Wiedemann type algorithms. So I can't help much with that.

If you want to do the alternating thing, there is an implementation here: https://github.com/thofma/Hecke.jl/blob/master/src/Sparse/UpperTriangular.jl
There is a function diagonal_form, which does the alternate thing and brings the matrix in diagonal form. Then one has to sort out the diagonal to get the Smith normal form.

By the way, do your HNF implementations also return the transform if asked for?

1. scalar_mul: nmod_mat and fq_mat_templates both have ordinary scalar mul with no suffix, with the former additionally having a scalar_mul_fmpz variant. fmpz_mat has fmpz, si, and ui variants, all of which are suffixed. I'm happy to have the sparse versions follow the convention of the corresponding non-sparse one, or always have a suffix, as you prefer.

One day we should make this consistent. There used to be a convention, which was confusing, so no one followed it. Probably just always having a suffix is the easiest convention to follow. You could open a ticket for this.

2. both fmpz_mat and fq_mat_templates define scalar_addmul and scalar_submul, but the operation is necessarily in place for fmpzs and either in or out of place for fq_mats. nmod_mat defines neither, but instead scalar_mul_add, which I think is the same as an addmul.

Are you referring to addmul/submul. I don't find the scalar ones for fq_mat. The scalar_mul_add should do a = ab + c not a = a + bc. If not, a ticket should be opened to rename this.

3. nmod_mat has a triple of functions get_entry, entry_ptr, and set_entry, which I replicated in the various sparse matrix libraries. However, fmpz_mat just has an entry function, with the same behavior as entry_ptr, while fq_mat_templates as both entry and entry_set, with the latter having the same behavior as set_entry.

Yes this is all inconsistent. You could open a ticket to make this more consistent. It's not urgent though.

4. set_nmod_mat in fmpz_mat has ordinary and unsigned modes, while the various CRT functions just take in a sign argument.

I'm not so worried about this one.

5. All three types of matrices have associated minpoly and charpoly functions, but for nmod_mat and fq_mat_templates, they are defined in the polynomial headers/folders, not in the mat headers/folders, which seems to be slightly unfortunate for namespacing. Is there a particular reason to put them there?

Probably there is, though there should be a note in the other header explaining why. If not, they can be safely moved, so long as everything compiles ok. Feel free to do this. I would refer they are in the mat header.

6. all three libraries have randrank and randops, which make an arbitrary diagonal matrix and do random operations on it respecitvely, but only nmod_mat combines them into a randfull.

Yes. Feel free to extend if you need it.

7. all three libraries have pretty printing, but fq_mat_templates and fmpz_mat also have regular printing (and file printing), and fmpz_mat has file/stdio reading.

Yeah, printing and file I/O are pretty hard. More functions are needed. Have at it if you are interested.

8. Random other omissions: fq_mat_templates is missing pow and trace, fmpz_mat is missing addmul, submul, and triangular solving (and random generation of triangular matrices), and neither nmod_mat nor fq_mat_templates have is_one, equal_row/col, gram, or kronecker product. 5, 6, 7, and 8 would be straightforward to "fix" if desirable. "Fixing" 2, 3, and 4 would break backward compatibility, though for (3) one could just introduce duplicate functions and keep everything compatible. And for 1, I'm happy to go either way as you like. K

Missing functions is just a measure of the amount of time developers have. It's much less serious than introducing arbitrary new conventions. Feel free to contribute to any of the above that interest you, and open (separate) tickets for the others.

I guess if you want to work on the breaking changes then it would be best to have those for the current release. If you could send me a list of name changes so that I can put them in the release announcement that would be helpful. Maybe we should add some #defines to the old names for now and deprecate them. Next release we can add some warnings to the deprecations and then the next release remove them.

1. scalar_mul: nmod_mat and fq_mat_templates both have ordinary scalar mul with no suffix, with the former additionally having a scalar_mul_fmpz variant. fmpz_mat has fmpz, si, and ui variants, all of which are suffixed. I'm happy to have the sparse versions follow the convention of the corresponding non-sparse one, or always have a suffix, as you prefer.

One day we should make this consistent. There used to be a convention, which was confusing, so no one followed it. Probably just always having a suffix is the easiest convention to follow. You could open a ticket for this.

Will do, will just follow the previous convention for the sparse mats for now.

1. both fmpz_mat and fq_mat_templates define scalar_addmul and scalar_submul, but the operation is necessarily in place for fmpzs and either in or out of place for fq_mats. nmod_mat defines neither, but instead scalar_mul_add, which I think is the same as an addmul.

Are you referring to addmul/submul. I don't find the scalar ones for fq_mat. The scalar_mul_add should do a = a_b + c not a = a + b_c. If not, a ticket should be opened to rename this.

In fq_mat_templates, lines 292 and 298 have scalar_addmul and scalar_submul, with C = A +/- c B (c the last argument). Addmul and submul (lines 273 and 279) are D = C + A*B.

In fmpz_mat, the scalar_add/submul are lines 212-220, and appear to be all B += cA or B -= cA (c the last argument). Divexact, and 2exp have the same in-place syntax. I don't see matrix addmul or submul for fmpz_mats.

In nmod_mat, scalar_mul_add(dest, X, b, Y) (declared on line 159, defined in scalar_mul_add.c) does appear to set dest = X + bY, i.e., the same operation as in fq_mat_scalar_addmul, but with a different order of arguments. Addmul and submul are also out-of-place, D = C + AB.

I went back and forth, and now have no opinion on whether out-of-place or in-place are better: any change would be trivial to make.

1. nmod_mat has a triple of functions get_entry, entry_ptr, and set_entry, which I replicated in the various sparse matrix libraries. However, fmpz_mat just has an entry function, with the same behavior as entry_ptr, while fq_mat_templates as both entry and entry_set, with the latter having the same behavior as set_entry.

Yes this is all inconsistent. You could open a ticket to make this more consistent. It's not urgent though.

Ok, I'll just add extra functions for now if that's ok, and we can decide later if anything wants to be deprecated.

1. set_nmod_mat in fmpz_mat has ordinary and unsigned modes, while the various CRT functions just take in a sign argument.

I'm not so worried about this one.

Cool, I copied the same format for sparse matrices.

1. All three types of matrices have associated minpoly and charpoly functions, but for nmod_mat and fq_mat_templates, they are defined in the polynomial headers/folders, not in the mat headers/folders, which seems to be slightly unfortunate for namespacing. Is there a particular reason to put them there?

Probably there is, though there should be a note in the other header explaining why. If not, they can be safely moved, so long as everything compiles ok. Feel free to do this. I would refer they are in the mat header.

I'll try to move them, and if it works, I'll do so in a new pull request.

1. all three libraries have randrank and randops, which make an arbitrary diagonal matrix and do random operations on it respecitvely, but only nmod_mat combines them into a randfull.

Yes. Feel free to extend if you need it.

Great.

1. all three libraries have pretty printing, but fq_mat_templates and fmpz_mat also have regular printing (and file printing), and fmpz_mat has file/stdio reading.

Yeah, printing and file I/O are pretty hard. More functions are needed. Have at it if you are interested.

Great.

1. Random other omissions: fq_mat_templates is missing pow and trace, fmpz_mat is missing addmul, submul, and triangular solving (and random generation of triangular matrices), and neither nmod_mat nor fq_mat_templates have is_one, equal_row/col, gram, or kronecker product. 5, 6, 7, and 8 would be straightforward to "fix" if desirable. "Fixing" 2, 3, and 4 would break backward compatibility, though for (3) one could just introduce duplicate functions and keep everything compatible. And for 1, I'm happy to go either way as you like. K

Missing functions is just a measure of the amount of time developers have. It's much less serious than introducing arbitrary new conventions. Feel free to contribute to any of the above that interest you, and open (separate) tickets for the others.

I guess if you want to work on the breaking changes then it would be best to have those for the current release. If you could send me a list of name changes so that I can put them in the release announcement that would be helpful. Maybe we should add some #defines to the old names for now and deprecate them. Next release we can add some warnings to the deprecations and then the next release remove them.

The only name changes would be the scalar_addmul/submul things, depending on what you want me to do: everything else is just adding new functions for now. Or does moving functions from one module to another also break things?

1. both fmpz_mat and fq_mat_templates define scalar_addmul and scalar_submul, but the operation is necessarily in place for fmpzs and either in or out of place for fq_mats. nmod_mat defines neither, but instead scalar_mul_add, which I think is the same as an addmul.

Are you referring to addmul/submul. I don't find the scalar ones for fq_mat. The scalar_mul_add should do a = a_b + c not a = a + b_c. If not, a ticket should be opened to rename this.

In fq_mat_templates, lines 292 and 298 have scalar_addmul and scalar_submul, with C = A +/- c B (c the last argument). Addmul and submul (lines 273 and 279) are D = C + A*B.

In which file? I can't see these anywhere.

In fmpz_mat, the scalar_add/submul are lines 212-220, and appear to be all B += cA or B -= cA (c the last argument).

Yes, that is correct. This is what addmul/submul should do.

Divexact, and 2exp have the same in-place syntax. I don't see matrix addmul or submul for fmpz_mats.

ok.

In nmod_mat, scalar_mul_add(dest, X, b, Y) (declared on line 159, defined in scalar_mul_add.c) does appear to set dest = X + b_Y, i.e., the same operation as in fq_mat_scalar_addmul, but with a different order of arguments. Addmul and submul are also out-of-place, D = C + A_B.

Ok, it looks like it was added by someone along with some other code. It's clearly incorrect.

Let's correct the order of arguments and make it scalar_addmul. Probably no one uses it anyway.

I went back and forth, and now have no opinion on whether out-of-place or in-place are better: any change would be trivial to make.

If it's possible to implement something in-place without allocating a new matrix, that can be a big saving for small problems. But it's only worthwhile if the algorithm permits it.

We really ought to have a convention like _inplace on the function name. But it's a bit late to change this.

In other modules in Flint (the polynomials for example) we implement an in-place function, then have the more general function call the in-place version without an extra allocation if there is aliasing of an input and an output.

But we can't do that for matrices, since in Flint the general matrix interface says that you should always provide a new matrix for the output to go in, unless otherwise documented. So basically people have to look up whether a function is in-place or not.

At the very least we should make sure that for example, lu is consistent across all dense modules though, and similarly consistent across all sparse modules.

1. All three types of matrices have associated minpoly and charpoly functions, but for nmod_mat and fq_mat_templates, they are defined in the polynomial headers/folders, not in the mat headers/folders, which seems to be slightly unfortunate for namespacing. Is there a particular reason to put them there?

Probably there is, though there should be a note in the other header explaining why. If not, they can be safely moved, so long as everything compiles ok. Feel free to do this. I would refer they are in the mat header.

I'll try to move them, and if it works, I'll do so in a new pull request.

Ok.

1. Random other omissions: fq_mat_templates is missing pow and trace, fmpz_mat is missing addmul, submul, and triangular solving (and random generation of triangular matrices), and neither nmod_mat nor fq_mat_templates have is_one, equal_row/col, gram, or kronecker product. 5, 6, 7, and 8 would be straightforward to "fix" if desirable. "Fixing" 2, 3, and 4 would break backward compatibility, though for (3) one could just introduce duplicate functions and keep everything compatible. And for 1, I'm happy to go either way as you like. K

Missing functions is just a measure of the amount of time developers have. It's much less serious than introducing arbitrary new conventions. Feel free to contribute to any of the above that interest you, and open (separate) tickets for the others.
I guess if you want to work on the breaking changes then it would be best to have those for the current release. If you could send me a list of name changes so that I can put them in the release announcement that would be helpful. Maybe we should add some #defines to the old names for now and deprecate them. Next release we can add some warnings to the deprecations and then the next release remove them.

The only name changes would be the scalar_addmul/submul things, depending on what you want me to do: everything else is just adding new functions for now. Or does moving functions from one module to another also break things?

Moving functions won't break anything, as far as I know. I mean, technically it could break someone's code. But I'm pretty sure they'll figure it out soon enough. If you like, add the changes to the ticket on deprecations (see below).

The scalar_mul_add definitely should be fixed. I have no idea how that slipped past. But send me a list of any other scalar functions you change, or better yet put them on the (recent) ticket which I think is called something like "Announce deprecations".

The current plan is for the first beta next Wednesday. But I will hold off until you are done with this one, as I think it is important enough.

1. both fmpz_mat and fq_mat_templates define scalar_addmul and scalar_submul, but the operation is necessarily in place for fmpzs and either in or out of place for fq_mats. nmod_mat defines neither, but instead scalar_mul_add, which I think is the same as an addmul.

Are you referring to addmul/submul. I don't find the scalar ones for fq_mat. The scalar_mul_add should do a = a_b + c not a = a + b_c. If not, a ticket should be opened to rename this.

In fq_mat_templates, lines 292 and 298 have scalar_addmul and scalar_submul, with C = A +/- c B (c the last argument). Addmul and submul (lines 273 and 279) are D = C + A*B.

In which file? I can't see these anywhere.

You're totally right, the header file has a placeholder section for scalar functions (line 273 of the header), so I must have filled them in when I did my first pass and forgotten about it. I also added in the addmul function, there was only a submul before, though it was out-of-place. So lets go in-place then everywhere, i.e. scalar_addmul/submul(B, A, c) performs B += cA or B -= cA, while addmul/submul (C, A, B) perform C += AB or C -= AB. If that sounds good to you, I'll take care of it in a new fork and make a pull request.

In fmpz_mat, the scalar_add/submul are lines 212-220, and appear to be all B += cA or B -= cA (c the last argument).

Yes, that is correct. This is what addmul/submul should do.

Divexact, and 2exp have the same in-place syntax. I don't see matrix addmul or submul for fmpz_mats.

ok.

In nmod_mat, scalar_mul_add(dest, X, b, Y) (declared on line 159, defined in scalar_mul_add.c) does appear to set dest = X + b_Y, i.e., the same operation as in fq_mat_scalar_addmul, but with a different order of arguments. Addmul and submul are also out-of-place, D = C + A_B.

Ok, it looks like it was added by someone along with some other code. It's clearly incorrect.

Let's correct the order of arguments and make it scalar_addmul. Probably no one uses it anyway.

I went back and forth, and now have no opinion on whether out-of-place or in-place are better: any change would be trivial to make.

If it's possible to implement something in-place without allocating a new matrix, that can be a big saving for small problems. But it's only worthwhile if the algorithm permits it.

We really ought to have a convention like _inplace on the function name. But it's a bit late to change this.

In other modules in Flint (the polynomials for example) we implement an in-place function, then have the more general function call the in-place version without an extra allocation if there is aliasing of an input and an output.

But we can't do that for matrices, since in Flint the general matrix interface says that you should always provide a new matrix for the output to go in, unless otherwise documented. So basically people have to look up whether a function is in-place or not.

At the very least we should make sure that for example, lu is consistent across all dense modules though, and similarly consistent across all sparse modules.

1. All three types of matrices have associated minpoly and charpoly functions, but for nmod_mat and fq_mat_templates, they are defined in the polynomial headers/folders, not in the mat headers/folders, which seems to be slightly unfortunate for namespacing. Is there a particular reason to put them there?

Probably there is, though there should be a note in the other header explaining why. If not, they can be safely moved, so long as everything compiles ok. Feel free to do this. I would refer they are in the mat header.

I'll try to move them, and if it works, I'll do so in a new pull request.

Ok.

1. Random other omissions: fq_mat_templates is missing pow and trace, fmpz_mat is missing addmul, submul, and triangular solving (and random generation of triangular matrices), and neither nmod_mat nor fq_mat_templates have is_one, equal_row/col, gram, or kronecker product. 5, 6, 7, and 8 would be straightforward to "fix" if desirable. "Fixing" 2, 3, and 4 would break backward compatibility, though for (3) one could just introduce duplicate functions and keep everything compatible. And for 1, I'm happy to go either way as you like. K

Missing functions is just a measure of the amount of time developers have. It's much less serious than introducing arbitrary new conventions. Feel free to contribute to any of the above that interest you, and open (separate) tickets for the others.
I guess if you want to work on the breaking changes then it would be best to have those for the current release. If you could send me a list of name changes so that I can put them in the release announcement that would be helpful. Maybe we should add some #defines to the old names for now and deprecate them. Next release we can add some warnings to the deprecations and then the next release remove them.

The only name changes would be the scalar_addmul/submul things, depending on what you want me to do: everything else is just adding new functions for now. Or does moving functions from one module to another also break things?

Moving functions won't break anything, as far as I know. I mean, technically it could break someone's code. But I'm pretty sure they'll figure it out soon enough. If you like, add the changes to the ticket on deprecations (see below).

The scalar_mul_add definitely should be fixed. I have no idea how that slipped past. But send me a list of any other scalar functions you change, or better yet put them on the (recent) ticket which I think is called something like "Announce deprecations".

The current plan is for the first beta next Wednesday. But I will hold off until you are done with this one, as I think it is important enough.

Yep, fine by me.

Now that the release is done, I can work on getting this merged. In the mean time, could you take a look at the test failures on appveyor. It seems fmpz_mat_mul_vec is defined twice, presumably once in mpoly and once in your code.

So it is declared in mpoly.h but only defined in fmpz_mat.h, which is why it doesn't break my make. Is that invalid under some gcc's? The syntax is identical:

fmpz_mat_mul_vec(fmpz *y, const fmpz_mat_t A, fmpz *x)

I think it is the Windows compiler that is complaining about it.

I think I will merge this branch into trunk once the tests are passing again.

The main thing I wanted to do from here is work out whether there are any better algorithms for the operations you support. I think we will eventually want to switch to dense algorithms for some of the problems when things become dense enough. You probably have a better idea than I do where that should happen.

As for the code and interface, I am happy with it. Better we get this into trunk and then improve it later if needed.

As for fmpz_mat_mul_vec it should be defined in the module it is declared, which looks like it should be in fmpz_mat given the name. People won't be able to find if it appears in multiple places or a module with a different name to the function.

The Windows compiler won't accept it in more than one .h file because it has to do something fancy to export symbols from dll's (that is what the FLINT_DLL is all about).

Please fix the copyright headers. I see a bunch of new files that erroneously claim to have been written by me 10 years ago :-)

Ah I am sorry about fmpz_mat_mul_vec. It should probably say something about aliasing and have the signature

void fmpz_mat_mul_vec(fmpz * v, const fmpz_mat_t M, const fmpz * u)

and it should be moved from mpoly/compose_mat.c to fmpz_mat/

Sometimes I put stuff in mpoly that shouldn't really be there.

@kartikv do you need help getting the tests to pass? Please feel free to put fmpz_mat_mul_vec in the proper place.

It looks like it can't find fmpz_sparse_mat_howell_form_mod referenced now. Perhaps this can be removed? Or is there some other issue.

I reconfigured, rebuilt, and retested everything in my setup and it works fine, so I'm not sure how to fix this (or what is wrong).

Ah I think the new modules just need to be added to CMakeLists.txt. That ought to fix it.

Sorry for the long delay, it's been a long year...I believe I have rebased from trunk in wbhart's repo, added the resizing functions (but no tests/documentation yet), and added myself as an author to those files. I'll deal with any bugs found by appveyor and add the extra tests/docs this week.

Fanstastic. Looking forward to some more testing/docs.

Welcome back - looking forward to having this merged.

@fredrik-johansson, do you think we can partially rescue this PR?

Kartik will be at the workshop if all goes well this time, presumably working on this.

Generics would be nice here, no?

certainly!

The matrices which were causing me problems (as report on the mailing list recently) are very sparse and so are good candidates for speeding up if a sparse implementation was available. So I hope so!